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ABSTRACT 

A statistical-stochastic model of the complete life cycle of North Atlantic (NA) tropical cyclones (TCs) is 
used to examine the relationship between climate and landfall rates along the North American Atlantic and 
Gulf Coasts. The model draws on archived data of TCs throughout the North Atlantic to estimate landfall rates 
at high geographic resolution as a function of the ENSO state and one of two different measures of sea surface 
temperature ( SST): 1 ) SST averaged over the NA subtropics and the hurricane season and 2) this SST relative to 
the seasonal global subtropical mean SST (termed relSST). Here, the authors focus on SST by holding ENSO to 
a neutral state. Jackknife uncertainty tests are employed to test the significance of SST and relSST landfall 
relationships. There are more TC and major hurricane landfalls overall in warm years than cold, using either 
SST or relSST, primarily due to a basinwide increase in the number of storms. The signal along the coast, 
however, is complex. Some regions have large and significant sensitivity (e.g., an approximate doubling of 
annual major hurricane landfall probability on Texas from —2 to +2 standard deviations in relSST), while other 
regions have no significant sensitivity (e.g., the U.S. mid-Atlantic and Northeast coasts). This geographic 
structure is due to both shifts in the regions of primary TC genesis and shifts in TC propagation. 


1. Introduction 

Intense tropical cyclones (TCs) are among the most 
devastating natural phenomena. Estimating the poten- 
tial for economic damage is a topic of high public in- 
terest and cuts across meteorology, climatology, and 
economics (e.g., Mendelsohn et al. 2012; Peduzzi et al. 
2012). Landfall risk assessments are used by the in- 
surance industry for setting rates and by governments 
for establishing building regulations and planning 
emergency procedures. Given the large coastal pop- 
ulations in harm’s way, estimating the short- and long- 
term evolution of TC hazard is crucial. 

There has been much interest in ascertaining and 
understanding trends in North Atlantic TC frequency 
and intensity and whether any such trends are related to 
anthropogenic climate change and/or natural climate 
cycles. North Atlantic (NA) sea surface temperature 
(SST) has risen in recent decades, driven in large part by 
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greenhouse gas forcing (Santer et al. 2006). Some sta- 
tistical studies have suggested a link between SST and 
NA TC activity (Emanuel 2005; Eisner et al. 2008), 
while other studies point to the spatial structure of SST 
as being more important than absolute SST (Vecchi and 
Soden 2007; Vecchi et al. 2008; Villarini et al. 2010; 
Ramsay and Sobel 2011; Villarini et al. 2012). Direct- 
numerical and downscaled-numerical transient climate 
simulations of the twenty-first century using the Intergov- 
ernmental Panel on Climate Change (IPCC) models 
produce widely varying TC results (Emanuel et al. 2008; 
Camargo 2013; Villarini and Vecchi 2012). The consensus 
among the climate models that are best at representing 
TCs is for reduced TC frequency by 2100 but an increased 
fraction of TCs in the highest intensity categories (cat; 
Bender et al. 2010; Villarini and Vecchi 2013). 

Whatever are the long-term trends in large-scale 
measures of TC activity, and however they are driven, it 
is ultimately necessary for risk assessment to estimate 
landfall rates at local geographic scales. This task is 
beyond the capability of global climate models. The 
most straightforward way to analyze landfall risk and its 
climate dependence is to use historical landfall events. 
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This approach is sound if there are a sufficient number of 
events, such as are found over large sections of coast and 
over many years. However, if one aims to study geo- 
graphic landfall rates at high geographic resolution and 
in addition use subsets of data years based on climate 
state, then sampling error becomes a major issue. This 
can be alleviated by exploiting data from outside the 
immediate region of interest, for example by supple- 
menting the local landfall wind data with overocean data 
(e.g., Emanuel and Jagger 2010). Another alternative is 
basinwide statistical models of the entire life cycle of 
TCs to generate sets of synthetic TCs much larger than 
the historical set and compute landfall probabilities 
from the synthetic sets. Various such models have been 
developed by private and public sector researchers 
(James and Mason 2005; Vickery et al. 2006; Emanuel 
et al. 2006; Rumpf et al. 2007). 

An advantage of a basinwide model is that it utilizes 
historical track information over the full basin: orders of 
magnitude more data than just at landfall. In effect, 
misses, as well as hits, inform the landfall analysis. A 
disadvantage is that while reliable landfall data go back 
into the nineteenth century, full-basin data are less re- 
liable prior to the mid-twentieth century. In addition, 
the added complexity of a full-basin model compared to 
a landfall model increases the possibility of model bias. 
Still, for coastal regions of historically few or no land- 
falls, a basinwide model is beneficial, as seen in the 
analysis of Hall and Jewson (2008), who show in out-of- 
sample tests that the increased precision outweighs any 
loss of accuracy. A second advantage to simulating the 
TC full life cycle is that landfall changes can be de- 
composed into changes in various TC components, such 
as genesis, propagation, and intensity. Here, the model 
initially described by Hall and Jewson (2007) is ex- 
tended and used to explore the relationships between 
regional landfall rates and SST. 

Section 2 reviews the data on which the model is built, 
and section 3 describes the statistical-stochastic model, 
focusing on components that have not been documented 
elsewhere. The model is evaluated in section 4, and the 
landfall-SST relationships explored in section 5. 

2. Data sources 

We use the Atlantic hurricane database (HURD AT) 
best-track data from 1950 to 2008 (Javinen et al. 1984) to 
build the statistical model. These data include position, 
central pressures, and maximum sustained wind speeds 
(k'max) every 6h along a TC’s trajectory. We restrict 
attention to named TCs. The date 1950 represents 
roughly the start of the era of routine aircraft re- 
connaissance. Moreover, Vecchi and Knudson (2011) 


have estimated that prior to 1950 the number of TCs 
missed by undersampling increases substantially. Im- 
portantly, all HURD AT storms that survive less than 
2 days are excluded from the analysis. There is evidence 
that changes in the frequency of such short-lived storms 
are spurious, due to underreporting earlier in the record 
(Villarini et al. 2011), and we wish to avoid such ambi- 
guity in the results. As part of the track propagation 
component of the model we use a composite annual 
cycle of 500-mb zonal wind, which has been constructed 
from 1950 to 2008 6-hourly National Centers for Envi- 
ronmental Prediction (NCEP) wind data. 

We employ SST and the Southern Oscillation index 
(SOI) as independent predictor variables, both of which 
are known to influence TC activity. In this study, we 
focus on SST and hold SOI to a neutral value for all 
experiments. Other factors have been shown to in- 
fluence NA TC activity as well [e.g., the solar cycle and 
the North Atlantic Oscillation (NAO); Eisner et al. 
2006; Eisner and Jagger 2008]. 

There has been considerable debate as to how SST 
should be viewed in relation to NA TC activity (Vecchi 
et al. 2008). The question is of enormous importance for 
long-range projections of TC activity, given the robust 
projections of increased SST over the next several de- 
cades. Recent work has emphasized the importance of 
spatial structure in SST, rather than SST directly, as 
being the best indicator for TC activity based on physical 
arguments and numerical simulations (Vecchi et al. 
2008). In particular, relative SST (relSST), the NA SST 
compared to the global tropical mean SST, has been 
shown to be at least as good a statistical predictor of NA 
TC counts as SST directly (Villarini et al. 2010; Vecchi 
and Knudson 2011; Villarini et al. 2012). There is no 
reason to expect a robust signal in relSST due to the ra- 
diative forcing of well-mixed greenhouse gases, though 
aerosol forcing is known to have affected the NA dis- 
proportionately in the twentieth century (Villarini and 
Vecchi 2012). One physical argument for relSST is that it 
is a proxy for potential intensity (PI; Vecchi and Soden 
2007): a measure of the upper limit of TC intensity 
(Holland 1997; Emanuel 2000) and a factor in TC genesis 
(Camargo et al. 2007). 

For an analysis of recent trends and near-term pro- 
jections, however, the choice between SST and relSST is 
less important than it is for long-range projections. In 
the near term, natural variability dominates anthropo- 
genic trends, NA SST and relSST are well correlated, 
and either one may be acting as a convenient proxy 
variable for modes of climate variability that affect TCs 
via several mechanisms. In this study, we perform sep- 
arate analyses using as independent variables SST di- 
rectly and relSST. Gridded SST data are obtained from 
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Fig. 1. Time series of NA MDR JASO mean SST (red) and NA 
MDR JASO mean relative SST (blue). 


the Met Office Hadley Centre (Rayner et al. 2003) and 
are averaged over the North Atlantic main development 
region (MDR) and the hurricane season July-October 
(JASO) to obtain annual values. For relSST, JASO 20°S 
to 20°N global averages are subtracted from MDR-JASO 
SST. As will be shown from the perspective of modeling 
current climate TC landfall rates statistically, use of SST 
and relSST yield similar results because SST and relSST 
are well correlated (Fig. 1). At several key points in the 
paper we show results for both SST and relSST. 

3. Methodology 

The goal of this study is to estimate the sensitivity of 
local landfall rates along the North American coast to 
NA SST and relSST. This is difficult using solely TC data 
at landfall because there are too few such data on a local 
basis. Instead we make use of an SST-sensitive statistical- 
stochastic model for TCs over the full NA (i.e., a sta- 
tistical model of TCs from birth to death that draws on 


basinwide data). In this way, much more data are ex- 
ploited than just landfall. In effect, information from 
data over the full basin is projected onto the coast. Once 
constructed, the model is run to generate “event sets” of 
synthetic TCs that are much larger than the historical set, 
and these synthetic TCs are used to compute landfall 
rates. Large event sets and corresponding landfall rates 
are generated for a range of fixed SST and relSST. The 
model components also depend on SOI; in this study we 
focus on the SST effects and hold SOI at a constant 
neutral value (zero anomaly value) for all experiments. 

The four components of the TC model are 1) genesis, 
2) propagation, 3) lysis, and 4) intensity. The genesis, 
propagation, and lysis components are similar to those 
described by Hall and Jewson (2007), with updates de- 
scribed by Yonekura and Hall (2011, 2013, manuscript 
submitted to J. Appl. Meteor. ClimatoL). The intensity 
component is new. The basic scheme for generating 
a synthetic TC is to first simulate a genesis event. Then, 
from the genesis location the TC is propagated in 6-hourly 
increments. At each 6-h position, the model determines 
whether or the not the storm suffers lysis (terminates), 
continuing until lysis occurs. For the TC intensity, a time 
series of the maximum sustained wind speed (F max ) is 
placed on the track. The model components are discussed 
individually below, with the most detail reserved for the 
new intensity scheme. 

a. Genesis 

The genesis component is identical to that described in 
Yonekura and Hall (2011) for the northwestern Pacific, 
which in turn is a generalization of the Hall and Jewson 
(2007) steady-state kernel density genesis to the case of 
climate-dependent variation. We perform local Poisson 
regression of the annual number of TC genesis events on 
SST or relSST and SOI. In a data circle of radius L cen- 
tered on location r. an amiual time series of genesis counts 
is Poisson regressed on SST and SOI, and the resulting 
Poisson rate is divided by the area of the circle to obtain 



Fig. 2. Genesis rates and mean annual count (inset numbers) for SST = (left) —2, (middle) 0, and (right) +2 in units of standard 
deviations in JASO MDR SST over 1950-2009. Contour intervals are 0.003, 0.006, 0.009, 0.012, 0.015, and 0.018 annual counts per 1-degree 
longitude-latitude box from light blue to dark red. White regions correspond to values below 0.003. 
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Fig. 3. As in Fig. 2, but for those dependent on relative SST. 


a rate per area. We then multiply by the area of a 1° box 
about r to obtain a grid of Poisson rate coefficients spanning 
the domain of interest (5°-55°N, 100°W-0°). During simu- 
lations the appropriate value of SST and SOI are multiplied 
by the local coefficients to obtain local rates. The Poisson 
distribution defined by this rate is then sampled randomly to 


obtain a count of genesis events for the grid box for the year. 
The counts are almost always 0 or 1, given the small region, 
though in principle higher counts can occur. 

The radius L of the local data circle in which historical 
counts are included for regression is determined by out-of- 
sample likelihood maximization and is equal to 700 km for 


105 115 125 135 145 




Fig. 4. The annual cycle of formation rates from an optimized space-date kernel density. Spatial fields are shown for 10-day windows 
centered on the day of year shown above each panel. For plotting, the kernel density is scaled to have the unit maximum, and the contour 
intervals are 0.1, 0.2, ... , 0.9, blue to red. White regions correspond to values below 0.1. 
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Fig. 5. Illustration of the dependence of track propagation on the annual cycle of 500-mb zonal winds, (a) His- 
torical tracks (blue) and a sample local regression region (red) in the Gulf, (b) Scatterplot of 6-houly zonal propa- 
gation speed vs the concurrent 500-mb annual cycle wind speed and the regression line (red), (c) Time series of the 
annual cycle zonal wind over the region, (d) Time series of the track speed predicted by the regression (red) and the 
actual track speeds (symbols). 


the 1950-2008 period. This course scale blurs the sharp 
latitudinal gradient of the formation rate at the southern 
flank of the MDR. The Poisson coefficient fields are “re- 
focused” by forcing the local rates to conform to a steady- 
state kernel density field calculated from the 1950-2008 
data, which has an optimized bandwidth of 210 km. In ef- 
fect, though the climate dependence of the genesis is less 
well resolved, the mean rates at least are constrained by 
the better-resolved historical averages. Genesis is pro- 
hibited over land on physical grounds. Overland power 
that results from the isotropic-averaging kernel is trun- 
cated and redistributed nearby over ocean. 

Figure 2 shows the genesis rate distribution for SST 
standardized anomalies of —2, 0, and +2. There is an 
approximate doubling of the basinwide mean formation 
rate from SST = —2 to +2 (7.6 to 14.6TCsyr _1 ). The rate 
does not increase uniformly. Rather, there is a greater 
increase in the MDR, resulting in an effective shift 
southeastward of the center of genesis activity. This 
genesis site shift has also been identified by Kossin and 
Vimont (2007) and Wang et al. (2011), who argue that it 


is a consequence of an eastward shift of the Atlantic warm 
pool [and by extension and increased MDR July- 
September (JAS) SST], as well as associated reductions in 
wind shear, possibly connected to the Atlantic meridional 
mode (AMM; Kossin and Vimont 2007). Figure 3 shows 
the same genesis maps but now with relSST as the in- 
dependent variable. The resulting sensitivity is similar to 
SST, with a 10% reduced overall formation rate change 
(8.1 to 14.2'rCsyr '). (There is a small component of 
the relSST series that is linearly related to SOI: SOI 
explains about 20% of the variance in relSST. When we 
use a “corrected” relSST signal — relSST with this linear 
component removed — the results are not significantly 
different than Fig. 3.) 

We also include seasonality in the genesis model (i.e., 
the date during the year the event occurred). This is ac- 
complished with a 3D (longitude, latitude, and date) ker- 
nel density. The date dimension of the kernel density is 
randomly sampled at each location to provide a date for 
the genesis event. Figure 4 shows cuts of the spatial density 
at 25 10-day windows spanning the hurricane season. 
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Fig. 6. (top left) The 1000 synthetic tracks for +2 (red) and —2 (blue) SST emanating from a single point, (top 
right) The SST +2 and —2 mean tracks (red and blue curves, respectively and the track-point density above 
a threshold (orange and light blue, respectively) about the +2 and —2 means, (bottom left) The +2 and —2 mean 
tracks (red and blue, respectively) for a 100-member jackknife uncertainty set. The means are significantly dif- 
ferent to the extent the distributions of means across the jackknife set are distinct, which they clearly are. (bottom 
right) Red (green) indicates regions of significant westerly (easterly) shift in mean track vectors with increased 
SST. In other regions, the propagation-SST relationship is not significant. Significance is determined by the same 
jackknife procedure. 


b. Propagation 

TC propagation is identical to that described by 
Yonekura and Hall (2011) for the NWP, which is a gen- 
eralization of the Hall and Jewson (2007) track model. 
Given a TC at location the two vector components of 
the 6-hourly increment dr to the new position r + dr are 
determined by linear multiple regression of nearby (to be 
determined objectively) HURD AT increments on SST 
(or relSST), SOI, and 500-mb zonal wind speed from an 
NCEP-based composite annual cycle. Regression on this 
composite, illustrated in Fig. 5, provides seasonal de- 
pendence of the propagation. The residuals are stan- 
dardized, assuming a generalized anisotropic correlated 
normal process, as illustrated in Fig. 3 of Yonekura and 
Hall (2011). The standardized residuals are modeled by 
lag-one autoregression [AR(l)]. HURDAT contribu- 
tions to the regression at r are weighted inversely with 
distance from r using a Gaussian kernel whose length 
scale is optimized by an out-of-sample likelihood max- 
imization. The 500-mb wind-independent variable is 
used because midlevel winds are known to steer TCs. 


The annual cycle provides a mechanism to reproduce 
seasonality in tracks. 

There is a statistically significant relationship between 
SST and track propagation in parts of the domain, as 
illustrated in Fig. 6. Tracks have large stochastic com- 
ponents, and launching many tracks from a fixed point in 
the middle MDR leads to a wide and overlapping spread 
of tracks, from which it is difficult to discern any impact 
of SST. When track-point density is contoured on the 
domain, there appears to be a slight shift away from the 
North American coast from cold to warm years. This 
shift is seen more clearly when just the mean tracks are 
plotted (i.e., the stochastic component is set to zero). 
However, alone it is not clear if this difference is sig- 
nificant. To test the significance we perform a general- 
ized jackknife test. The regression of HURDAT tracks 
on the independent variables is performed 100 times, 
each time dropping a random 14 (—20%) of the 69 data 
years. (The fraction 20% is an arbitrary compromise 
between two competing factors: the fraction should be 
large so that there are many possible permutations for 
the subset years from the full set, and the fraction should 
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Fig. 7. Root-mean-square forecast error in maximum V max along 
tracks as function of SST anomaly scale. A shallow minimum 
occurs at S = 0.5. 


be small so that there are many data within the subset to 
perform regressions.) For each 55-yr subset, the mean 
tracks are simulated from the selected point for SST = +2 
and —2. The result is shown in Fig. 6 (lower left). 
Clearly, the jackknife hot and cold sets are distinct, and 
the shift in mean track with SST is significant, at least 
from this selected location. Finally, the jackknife pro- 
cedure is repeated for each location on a 1° by 1° grid, 
using single-step (6 h) mean tracks. The lower right 
panel of Fig. 6 shows regions where the SST depen- 
dence is significant: red indicating an eastward shift in 
propagation and green indicating a westward shift in 
propagation. 

Wang et al. (2011) find a similar effect of SST on track 
propagation using dynamical simulations. They argue 
that a large Atlantic warm pool (and by extension 
a warmer MDR JAS mean SST) leads to a stationary 
wave pattern that weakens the NA subtropical high. A 
weakened subtropical high in turn induces an eastward- 
steering flow anomaly in the western Atlantic off the 
North American coast. Similarly, Colbert and Soden 
(2012) note an increase in recurving versus straight- 
moving TCs from the MDR with a weakened sub- 
tropical high, a change in propagation that they also 
associate with the AMM. The weaker relationship we 
find between track and relSST (not shown) may simply 
be a signature of the track-SST relationship, via the 
SST-relSST correlation. 

c. Lysis 

At each 6-h position r, along a simulated track, the TC 
has a probability of terminating. The probability P is 
derived from analysis of the historical TCs that have 


terminated near r. as described in Flail and Jewson 
(2007): 

I e~^ IL 


all 

where the sum in the denominator is taken over all TC 
points, and the sum in the numerator is taken only over 
terminal TC points. The length scale L determines how 
locally the sums should be weighted and is determined 
by out-of-sample log-likelihood maximization. The 
meteorological mechanisms for TC dissipation over 
water and land are very different. Therefore, land and 
ocean lysis probabilities are strictly separated using 
a 0.1° land-ocean mask. 

Note that in formulation (1) the TC’s intensity does 
not appear. Flowever, the time series of F max that gets 
placed on simulated tracks (see below) is constrained by 
the lysis location such that intensity at lysis is realistic. 

d. Intensity 

The model component for time series of the maximum 
sustained wind speed (F max ) is different in character 
than the other components. Instead of regression, we 
employ a scheme to resample the HURD AT F max time 
series. First, a track is simulated using the genesis, 
propagation, and lysis schemes described above. Then, 
a U max series is selected among the HURD AT series by 
a weighted random draw, and this series is placed on the 
simulated track, with rescaling in time to match the 
duration of the simulated track. Finally, a random per- 
turbation to the maximum K max along the series is made, 
and the entire series is rescaled to match the new per- 
turbed maximum F max but holding the first and last 
points unchanged. 

The random draw of a HURD AT TC’s F max series is 
weighted toward HURDAT TCs whose tracks are 
similar to the simulated track on which the F max series is 
to be placed. Similarity is defined by several criteria: 1) 
the proximity of the simulated and historical genesis 
sites, 2) the proximity of the simulated and historical 
lysis sites, and 3) the proximity of the amount of time 
that simulated and historical TCs spend in different 
latitude bands. If a track makes one or more landfalls, 
then it is broken into overocean and overland segments, 
and separate weighted U max random draws are made for 
each segment. In such landfalling cases additional sam- 
pling weights are used: 4) the proximity of the simulated 
and historical landfall locations and 5) the proximity of 
the last overocean F max value to the first overland F max 
value. If the simulated TC reemerges over the ocean, 


1 November 2013 


HALL AND YONEKURA 


8429 


(a) (b) 



Fig. 8. (a)-(c) Three randomly selected examples from 1950-2010 simulations, (d) along with the HURD AT tracks. Blue corresponds to 

cat-0-2, yellow cat-3-4, and red cat-5. 


then the new V max series are selected for the reemergent 
TC, again using proximity weights at the point of re- 
emergence. The transition from one V max series to an- 
other ocean to land and land to ocean is made smooth by 
rescaling the new segment to match the prior segment at 
the transition point. 

We want to allow for direct sensitivity to the inde- 
pendent climate variables in the V max selection. (There 
is already indirect sensitivity. The genesis sites change 
with SST and SOI, and different Vmax series will be 
preferentially selected.) Therefore, a final selection 
weight is added: the proximity of the SST (or relSST) 
anomaly from the year of the candidate V max series to 
the SST (or relSST) anomaly of the current simulation 
year. That is, in the random draw a candidate V max series 
i gets weighted by exp[— (SST s ; m — SST,)/L], where 
SST S i m is the SST anomaly of the simulation year, SST, is 
the SST anomaly of the year of the candidate V max se- 
ries, and L is a scale to be determined. (The quantity L is 


dimensionless, as the SSTs are expressed as standard- 
ized anomalies.) 

Does the addition of such an SST weight lead to im- 
provement in model forecast or overconstrain the model? 
To answer this, we perform an out-of-sample year opti- 
mization of the scale L used to define SST proximity. On 
the one hand, if it turns out that the optimal L is infinite, 
then V max series from all years, regardless of SST (or 
relSST), should be used equally, and preferential weight- 
ing by SST is overconstraining. On the other hand, if the 
optimal L is finite, then the SST (or relSST) state is im- 
portant enough to warrant selecting from a V max dataset 
effectively smaller than the full set. 

The procedure is as follows: 

1) Pick an SST anomaly scale L. 

2) Pick an out-of-sample year in 1950-2008 and loop 
through the historical tracks of that year. 

3) For each track, select a V max series from the remain- 
ing years. 
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Fig. 9. Map illustrating 87 coast segments along the Caribbean, Atlantic, and Gulf Coasts. The segments are approximately 100 km in 
length. They are divided into 10 regional categories, through color codes and labels. 


4) Compute the difference of the maximum V max of the 
selected series from the maximum V max of the 
historical series. 

5) Accumulate the squared differences and record the 
root-mean-square average as the out-of-sample fore- 
cast error for the anomaly scale. 

6) Minimize the error in L. 

This type of optimizing is similar to the optimization 
performed in the selection of scales defining the locality 
of regression in the track, genesis, and lysis models 
(Hall and Jewson 2007). There is a balance between 
two competing factors: 1) a large L includes as much 
data as possible, thereby reducing sampling error; and 
2) a small L resolves underlying meteorological de- 
pendence of TC intensity on SST. Figure 7 shows the 
error as a function of the SST anomaly scale. There is 
a shallow minimum at L = 0.5, indicating that V max has 
enough sensitivity to SST to warrant the inclusion of 
SST as a weight, despite the increased sampling error. 

In this resampling scheme for V max , the maximum 
intensity possible for a simulated storm is the maximum 


intensity recorded in the HURDAT database. This is 
a limitation and is not in keeping with the stochastic 
nature of the other model components, which are not 
bounded by historical values. To remedy this we add 
random perturbations to the maximum V max values 
achieved along TCs [lifetime maximum intensity 
(LMI)], rescaling each storm’s V max series accordingly, 
such that the first and final V max values are unchanged. 
This rescaling, however, needs to be done carefully, so 
that the distributional character of the resulting set of 
LMI is unchanged. We do not want the resulting set of 
simulated LMI to have a frequency distribution differ- 
ent in character than the historical set of LMI. 

First, we perform a simulation and store the set of as 
yet unperturbed LMI. This unperturbed simulation set 
has a frequency distribution of LMI that matches the 
historical distribution, assuming the track simulations 
are realistic, so that an unbiased sample of the V max 
series is selected. Following Casson and Coles (2000), we 
model the LMI with a generalized extreme value (GEV) 
distribution fit to the historical LMI. GEV is the ap- 
propriate distribution for the maximum in a series of 
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Fig. 10. Landfall return-period curves for the full coast and each of the 10 regions, as labeled. Yellow curves are shown for each of the 
5000 simulations of the 59-yr (1950-2008) period. The orange curves indicate the lower and upper bounds of the inner 95% across the 5000 
simulations. The red curve is obtained by combining the 5000 simulations in series (i.e., a 295 000-yr simulation). The blue curve is 
obtained directly from HURD AT landfalls. Model landfall is unbiased to the extent that the blue curve falls inside the orange curves. 


random values. The GEV distribution is then converted 
to normal by matching cumulative distribution functions 
(CDFs). (That is, each maximum V max value in the GEV 
has a CDF value Y. The quantity X is found, such that 
the standard normal CDF at X equals Y.) The elements 
of the normal set are perturbed with a normal random of 
zero mean and standard deviation cr, which is de- 
termined objectively by an out-of-sample minimization 


of the LMI forecast error. The perturbed normal set is 
then standardized. Finally, the standard normal per- 
turbed set is converted back to GEV, again by CDF 
matching. This results in a set of V max series on simu- 
lated tracks whose maximum V max is not bounded by 
the highest observed V max but whose frequency distri- 
bution has a form that matches that of the historical 
maximum F max . 
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category 

Fig. 11. Annual landfall rates (counts yU 1 ) over the full coast (Yucatan through U.S. Northeast) of Fig. 10 
separately in each intensity category 0 to 5. The red symbols are the mean of the 5000-member 1950-2008 
simulation set. The orange bands show the 95% ranges about these means. The blue symbols are the rates 
directly computed from HURD AT. 


4. Evaluation 

Before using the TC model to estimate landfall rate 
sensitivity to SST and relSST, it is important to evaluate 
the model. The 1950-2010 period is simulated 5000 
times and the results compared to the historical TCs. 
Tracks from three such simulations are shown in Fig. 8, 
along with the historical tracks. The model performs 
well according to a specified diagnostic to the extent that 
the historical diagnostic value appears as a typical 
member of the much larger synthetic set of diagnostic 
values. The focus here is hurricane landfall, and we use 
landfall rates to evaluate the model. For each simula- 
tion, landfalls are computed in 100-km segments along 
the Atlantic-Caribbean coast from the Yucatan Pen- 
insula through Maine (Fig. 9). For each 1950-2008 
synthetic set the landfalls are used to construct return- 
period curves, which are the average times between 
successive landfalls above a U max threshold as a function 
of the threshold. Because of the stochastic nature of the 
simulations, a return-period curve is a statistic, and no 


two return-period curves are identical. These synthetic 
curves can be compared to the return-period curves 
computed directly from historic landfalls over the same 
period. The model is unbiased if the historic return- 
period curve falls within the inner 95% of the simulated 
curves. Figure 10 shows these return-period curves, with 
one panel for each of the 10 regions as labeled. In most 
regions at most intensity thresholds, the model is in- 
deed unbiased by this standard. There is some evidence 
for a low bias at low V max thresholds in the northern 
Gulf Coast. 

Figure 11 shows the landfall rates (counts yU 1 ) along 
the full coast, Yucatan to the U.S. Northeast, in each 
intensity category separately. Shown are the model 
means across the 5000-member 1950-2008 simulation 
set, the 95% band about this mean, and the rates com- 
puted directly from HURDAT. For all categories but 
one, the HURDAT rate falls well inside the uncertainty 
range about the model mean, indicating the model 
landfall over the full coast is not biased for these in- 
tensities. For cat-1 landfalls, the HURDAT rate sits 
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Fig. 12. Probability forecasts of NA TC counts in 2009 (blue curve) and 2010 (red curve) for 
(a) all named TCs and (b) major (cat-3+) hurricanes. The symbols indicate the actual numbers 
that occurred in 2009 (blue) and 2010 (red). 


right at the top of the range, indicating a marginal low 
model bias for this category. 

We further evaluate the model with an out-of-sample 
test of its forecast ability. The model has been con- 
structed using HURD AT 1950-2008, but the 2-yr period 
2009/10 has also been included in the simulations. We 
compare these simulations to 2009/10 observations. In 
other words, given observed 2009 and 2010 values for 
the independent variables SST and SOI, we test the 
model’s probabilistic forecast of these years, which were 
not used in its construction. (A full forecast would re- 
quire predictions for SST and SOI, as well; here, just the 
TC model is being tested.) The year 2009 was moder- 
ately unfavorable for NA TCs, with an SST anomaly of 
—0.4 and an SOI anomaly of 0.9 (El Nino), while 2010 
was extremely favorable, with an SST anomaly of 3.1 
and an SOI anomaly of —1.9 (La Nina). 

Figure 12 shows the 2009 and 2010 distributions over 
the 5000-member simulation ensemble of NA-wide 
count of (i) all named TCs and (ii) all TCs that some- 
where achieve major hurricane status (cat-3+). 
Although there is overlap between the 2009/10 distri- 
butions, there is a clear shift to more TCs in 2010. The 
ensemble-mean total TC counts are 9.6 (2009) and 18.3 
(2010), and the ensemble-mean major hurricane counts 


are 2.1 (2009) and 5.8 (2010). The numbers that occurred 
in these years, shown as symbols in the figure, are 9 
(2009) and 19 (2010) for total TCs and 2 (2009) and 5 
(2010) for major hurricanes. These values are right at the 
peaks of the forecast probabilities, indicating good 
model forecasts. 

Figure 13 shows the 2009/10 distributions over the 
5000-member simulation ensemble of Yucatan-to-Maine 
landfalls of all named TCs and major hurricane status 
(cat-3+). The model does equally well for landfall fore- 
cast as for basinwide storm count. The mean forecast for 
all landfalls in 2009 is 3.8 and in 2010 is 7.2. The actual 
occurrence was 2 and 8, both values near the peaks of 
their probabilistic forecasts. The mean forecast for major 
hurricane landfalls in 2009 is 0.4 and in 2010 is 1.1, while 
the actual occurrence was 0 and 1 (Hurricane Karl in 2010 
on the Mexican Gulf Coast), again near the peak prob- 
abilistic forecast. 

5. Landfall and SST 

The primary goal of this study is to estimate the re- 
lationship between SST (or relSST) and landfall rates. 
For each of the five fixed values of SST anomaly (—2, — 1, 
0, +1, and +2) we perform 10000-yr simulations. The 
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Fig. 13. As in Fig. 12, but for landfalling TCs along the full segmented coast shown in Fig. 9. 
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Fig. 14. (top) Probabilities along the U.S. and Mexican coasts (as 
labeled) of one or more TC landfalls of any intensity (cat-0+) per 
year per 100 km of coastline in +2 sigma SST years (red) and 
-2 sigma SST years (blue), (bottom) Difference in such proba- 
bilities hot minus cold. In each curve, the spread is a measure of 
significance, representing the inner 95% of results in a generalized 
jackknife test. In this figure, all of the model components genesis, 
propagation, and intensity are sensitive to SST. 


long duration is necessary for the mean landfall rates at 
the 100-km resolution to converge sufficiently and long 
return period to be estimated. We then compare mean 
landfall rates along the coast between the extreme +2 
and —2 states. 

It is crucial to estimate the uncertainty of the mean 
rates. A landfall rate difference between two SST states 
may be found, but is it significant? Uncertainty bounds 
about the means are determined with a generalized 
jackknife test. For each value of SST we perform 100 
versions of the 10000-yr simulation. In each version, all 
model components are reconstructed on a random 80% 
subset of the data years 1950-2008 (i.e., a random 47 of 
the 59 years). This tests the sensitivity of the results to 
the finite nature of the underlying data on which the 
results are based. With this test we can say that the dif- 
ference in mean landfall rates on some coastal region 
between the SST = +2 and SST = —2 states are sig- 
nificantly different (or not) if the inner 95% of the set of 
100 differences excludes (or does not exclude) zero. 

Figure 14 shows the probability of TC landfall at any 
intensity per year per 100-km along the coast from the 
eastern Yucatan to Maine for warm (+2) and cold (—2) 
years. Also shown is the hot-cold difference. There is 
considerable geographic structure in the probabilities. 
In either hot or cold years, the eastern Yucatan and the 
northern Gulf Coasts have the greatest exposure, while 
the western Yucatan has almost no probability of a di- 
rect landfall. (Storms crossing the Yucatan east to west 


Fig. 15. As in Fig. 14, but for major hurricanes (cat-3+). 

can do heavy damage on the west side of the peninsula. 
The diagnostic here counts storms as they actually cross 
ocean to land.) Figure 15 shows the same landfall prob- 
abilities but now for major hurricane landfalls (cat-3 + ). 
The geographic structure is similar to TC landfall at any 
intensity, though there is increased weight toward the 
Gulf and away from the U.S. East Coast and Northeast. 
Figure 16 shows the same major landfall probabilities 
but now with relSST as the independent variable. The 
results are very similar to SST, the biggest difference 
being western Florida, where the warm-year increase is 
marginally significant for SST but is insignificant for 
relSST. 

With increased SST (or relSST), there is increased risk 
of landfall on the eastern Yucatan, on Texas, on the 



Fig. 16. As in Fig. 15, but the independent variable is relative 
SST (NA MDR JASO mean minus global northern subtropical 
JASO mean), with red +2 sigma and blue —2 sigma. 
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Fig. 17. Return-period curves for the fixed SST states +2 (red) and —2 (blue). The spread for each state indicates 
the 95% confidence bounds from the jackknife uncertainty test. The first panel includes the entire segmented coast of 
Fig. 9, while subsequent panels represent the smaller coastal regions of Fig. 9, as labeled. Note the different year 
range on the x axis for the FULL COAST panel. 


north Gulf Coast, and (for SST) marginally on the 
Florida Gulf Coast. The biggest increase is seen on the 
eastern Yucatan, where the probability of a cat-3 + 
landfall is 2 to 3 times larger for +2 SST than —2 SST. In 
Texas, the increase in landfall rate from —2 to +2 SST 
(or relSST) is roughly a factor 2. Interestingly, the U.S. 
East Coast from Florida through New England shows no 
significant sensitivity to SST or relSST. This is explored 
further below. 

Figure 17 shows the SST = +2 and —2 landfall-Umax 
return-period curves for the full coast of Fig. 9 and for 
the individual regions. The spread for each SST state 
corresponds to the jackknife 95% confidence bounds. 
The eastern Yucatan coast has a reduction in return 


periods (increase in landfall frequency) that is significant 
at all intensities. That is, the red and blue regions are 
distinct at all intensities. Other Gulf Coast regions dis- 
play a more varied SST sensitivity, typically showing 
a significant reduction in warm-year return periods for 
some intensities but not others. For example, Texas has 
significant return-period reduction for V max less than 
40 kt (lkt = 0.514 44 ms -1 ) and V max between 80 and 
120 kt. This complex signal is due to the tendency for 
warm-year return-period curves on the Gulf to become 
more vertical in the cat-2-3 intensity range, consistent 
with the observation that cat-2 and cat-3 landfall rates 
differ by relatively little (see cat-2 and cat-3 landfall 
rates in Fig. 11). The U.S. East Coast, from Florida 
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Fig. 18. As in Fig. 15, but only genesis is sensitive to SST. 

through New England, shows no significant changes in 
landfall return periods at any intensity. Summed up for 
the full coast, there are significant reductions in return 
periods for all intensities up to about 120 kt, except for 
a small range around cat-2. 

To help interpret the different sensitivities to SST 
along the coast, the entire analysis is repeated but now 
making only subsets of genesis, propagation, and in- 
tensity sensitive to SST. Figure 18 is similar to Fig. 15, 
except that only genesis is sensitive to SST. (That is, 
there is no SST regression in the propagation compo- 
nent, and the F max sampling scheme now considers all 
years equally, independent of SST.) With only genesis 
sensitive to SST, there is marginally significant landfall 
increase everywhere along the Caribbean and Gulf 
Coasts through Florida. The increase on the U.S. East 
Coast is marginally insignificant. If TC genesis increased 
uniformly over the NA with SST (or relSST) then, 
without SST dependence in any other model compo- 
nent, we would expect the landfall rate change to be 
uniform along the coast. The modest geographic struc- 
ture observed in the genesis-only landfall rate change 
must be because NA genesis does not change uniformly. 
Instead, the center of genesis shifts southeastward in 
warm years (Figs. 2 and 3). 

When SST sensitivity is added to the track propaga- 
tion (Fig. 19), the landfall profile looks more like that of 
the full model. The Yucatan shows no additional effect, 
the warm-year increase on the Mexican Gulf and Texas 
is exaggerated, and the rest of the Gulf and the U.S. East 
Coast no longer shows SST sensitivity. Figure 6 shows 
that in warm years TCs originating in the MDR do not 
reach as far west before recurving northeastward. This 
reduces the odds of a U.S. East Coast landfall, buffering 
the region from the overall increased TC count. Finally, 



Fig. 19. As in Fig. 15, but only genesis and tracks (not F max ) are 
sensitive to SST. 


with F max also sensitive to SST or relSST (Figs. 15 and 
16), the Yucatan and Mexican Gulf sensitivities are ex- 
aggerated, Texas and the north Gulf have significant 
landfall increase with SST or relSST, and western Flor- 
ida has significant increase with SST (but not relSST). 
The U.S. East Coast, however, remains insensitive to 
SST or relSST, mostly due to the buffering effect of track 
propagation changes, as well as the shift southeastward 
in the region of primary genesis. 

To further illustrate the impact on landfall of factors 
beyond changes in basinwide counts, we plot in Fig. 20 
the fraction of TCs of any intensity that make landfall 
per 100 km along the coast. That is, we divide the 



km 


Fig. 20. (top) The fraction of all NA TCs that make landfall per 
100 km of segmented coastline for hot (red) and cold (blue) SST. 
The spread indicates the inner 95% of the jackknife set. (bottom) 
The inner 95% of the hot-cold difference of the landfall fraction. 
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Fig. 21. Normalized histograms across the jackknife set of the hot-cold differences in landfall fraction of all NA TCs. Blue shows the 
inner 95% and red shows the upper and lower 2.5% across the jackknife set. (a) The entire coastline, Yucatan to the U.S. Northeast, 
(b) The Yucatan and Mexican Gulf, (c) Florida through the Northeast. 


landfall rates in the hot and cold SST simulations by the 
total hot and cold basinwide TC count. The landfall 
fraction depends on changes in genesis location and 
track with SST but not overall formation rate. From 
Florida through the Northeast there is a marginally 
significant decrease in landfall fraction with high SST 
compared to low. In contrast, from the Yucatan through 
the Mexican Gulf there is an increase warm over cold; 
the increase is significant on the eastern Yucatan coast 
but marginally insignificant on the Gulf Coast according 
to our jackknife uncertainty criterion. 

Figure 21 summarizes these results in terms of the 
normalized histograms over the jackknife set of the 
difference of hot over cold landfall fraction, now sum- 
med over (i) the entire coast, (ii) the Yucatan through 
Mexican Gulf, and (iii) Florida through the U.S. 
Northeast. For the entire coast, there appears to be 
a decrease in landfall fraction with high SST, but the 
difference is not significant. The apparent increased 
landfall fraction on the Yucatan and Mexican Gulf is 
also not significant. The decrease on the region Florida 
through the Northeast is significant. In summary, we find 
no significant change with SST in the fraction of NA TCs 
making landfall anywhere from the Yucatan to the U.S. 
Northeast, in agreement with Villarini et al. (2012). 
However, the combined effects of shifts in genesis site 
and track changes hint at a dipole landfall effect, with 
a marginally insignificant increase in warm years on the 
Yucatan and western Gulf Coast, and a marginally sig- 
nificant decrease on the U.S. East Coast. 

These results on the geographic distribution of land- 
fall rate changes with SST are at least qualitatively in 
agreement with the work of Kossin et al. (2010). These 
authors decomposed NA TC tracks into four clusters 
and regressed the annual formation rate of each cluster 
on the AMM, which has a strong signature in NA SST. 
The cluster of TCs originating in the MDR and propa- 
gating westward to make western Gulf and Yucatan 


landfall is particularly sensitive to AMM. Its formation 
rate increases with increasingly positive AMM (and 
hence NA SST), leading to the enhanced western Gulf 
landfall rates we see in our analysis. 

6. Conclusions and discussion 

We have employed a statistical-stochastic model of 
TCs in the North Atlantic to estimate the sensitivity of 
North American TC landfall to SST and relative SST 
(relSST). We find statistically significant increases in 
annual major hurricane (cat-3 + ) landfall probability in 
warm years compared to cold for the eastern Yucatan, 
the Mexican Gulf Coast, Texas, and the north Gulf 
Coast. In addition, there is a significant increase on 
western Florida for SST, though not for relSST. In- 
terestingly, there is no significant landfall change on the 
U.S. Atlantic coast, despite the overall increased annual 
TC count. Selective removal of the SST dependence 
from key model components reveals that the U.S. East 
Coast is buffered from the increased TC count by a shift 
with SST (and relSST) in track propagation from TCs 
originating in the MDR. In warm years, these tracks 
recurve to the northeast on average at a point farther 
east than in cold years, thereby reducing their odds of 
making East Coast landfall. Also contributing to this 
East Coast buffering is southeastward in the center of 
NA TC genesis and a consequent reduction in the 
landfall odds. 

A purely statistical study such as this can make no 
statement about the physical mechanism of the SST- 
landfall relationship. Kossin et al. (2010) show statisti- 
cally that the SST patterns favorable for MDR TCs are 
high SST in the NA and low SST in the subtropical Pa- 
cific. Motivated by this and other studies we also tested 
relSST, in place of SST, as the independent predictor 
variable. The landfall results change little. Over the past 
60 years, apparently, the climate state when SST is high. 
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either absolutely or relatively, is associated with NA TC 
changes that cause significant landfall increases on parts 
of the North American coast. SST and/or relSST may be 
acting as a proxy for other TC-influencing variables as- 
sociated with natural climate modes. For example, 
Kossin and Vimont (2007) discuss the Atlantic meridi- 
onal mode (AMM) and its impact on the Atlantic warm 
pool and wind shear. In the AMM’s positive phase, the 
North Atlantic is warm, providing a link to our SST and 
relSST predictors. The eastward expansion of the warm 
pool shifts the center of genesis eastward, which affects 
landfall rates, as TCs originating in the eastern MDR 
have lower odds of making North American landfall. 
Wang et al. (2011) have also argued that, via standing 
wave influence on the subtropical high-pressure system, 
the expanded warm pool influences TC propagation in 
precisely the way seen here: less westward penetration 
before northeastward recurvature, which reduces land- 
fall odds. 
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